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SUMMARY 


The methods of Hamilton and Lagrange in the calculus of 
variations are employed to optimize planetocentric orbital transfers 
using fixed-thrust engines, where the thrust acceleration ¢ ¢ local 
gravitational acceleration. The equations were derived for minimum 
time transfers. However, since fixed-thrust engines are used, the 
minimum fuel problem is solved simultaneously. Numerical solutions are 
obtained by integrating the optimized differential equations of motion 
for the problem, the actual integrations being performed by an IBM 360/30 
computer. The numerical integrations result in optimal spiral 
trajectories representing orbital transfers between specified coplanar, 
co-axial ellipses, with the optimal thrust angle being calculated at all 
points along the trajectories. For an orbital transfer at the planet 
Mars, the low fixed-thrust system yields a significant increase in the 
amount of payload that can be delivered into the desired final orbit 
compared with a transfer accomplished by an impulsive (high thrust) chemical 
rocket propulsion system. However, the time required to achieve the 
transfer is necessarily longer in the low-thrust case. It is then necessary 
to examine the trade-off between the increase in payload and the extra 
transfer time to determine whether or not the fixed-thrust system is 


desirable. 
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I. INTRODUCTION 

we we explore the space around us, we look for ways of doing things 
more efficiently, ways to save money, and ways to put heavier payloads into 
orbit at the same or less cost than is presently possible. One feasible way 
of accomplishing this is through the use of solar-electro-static or nuclear- 
electro-static propulsion systems. They require much less fuel (weight- 
wise) than chemical rocket engines, have much higher specific impulses, and 
can result in considerable increases in final payload delivered. They do, 
however, require more time to achieve the same orbital maneuvers because of 
the much lower thrust inherent in these type of systems. 

Considerable research on optimal low thrust maneuvers has been done 
in recent years. A large part of this work has been done by T. N. Edelbaum. 
He has published papers which contain solutions for optimum power-limited 
transfer and rendezvous between arbitrary coplanar and coaxial elliptic and 
circular orbits. The assumptions involved in the development are 1) the 
thrust acceleration is much Smaller than the local gravitational acceleration, 
2) the propulsion system is power-limited, i.e., it operates at constant 
exhaust power with a variable thrust magnitude inversely proportional to the 
exhaust velocity. 3) the direction and magnitude of the thrust, both of 
which are assumed to be completely variable, is to be determined as a function 
of time so as to minimize the fuel consumption for a fixed total transfer time. 
Edelbaum has also developed an analytic solution for fixed thrust systems, but 


this is only applicable to transfers between almost circular orbits, for small 
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changes in a and e , and for large changes in a and i (orbital plane 
inclination)(Reference 1, 2, 3). 

This thesis addresses the problem of optimizing fixed, low- 
thrust orbital transfers with large changes in a and e. Based upon a 
Lagrange-Hamilton formulation developed by Dr. M. Handelsman (Reference 4) 
the numerical solution of optimized transfer trajectories is performed by 
computer. The solutions involve the optimization of planetocentric transfers 
between coplanar, coaxial ellipses or circles of specified characteristics, 
with no time constraint. Because of the relatively small fuel expenditure 
in these transfers, it is possible to replace "fixed thrust'' with "fixed 
thrust acceleration" throughout the entire transfer with little error. Solar 
electric engines lose a certain amount of thrust as their distance from the 
Sun increases; this could hardly be considered constant thrust in heliocentric 
Space. But, since the problem has been limited to a planetocentric study, 
the thrust of the solar electric engine can be assumed to be constant within 
the planetocentric space. 

By fixed thrust, it is meant that the thrust magnitude is constant 
and that the engine remains turned on at all times during the transfer. The 
problem then remains, where to point the thrust at all times. Now, there are 
many factors which could be optimized in an orbital transfer of this type. 

In this thesis two separate phenomena were optimized. First, orbital trans- 
fers were optimized for minimum time using fixed thrust engines. Since the 
thrust is fixed, and, therefore, continuous, then a minimum time transfer 
would also be a minimum fuel transfer, the amount of fuel expended being 


directly proportional to the length of time the engine is turned on. Secondly, 
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orbital transfers were optimized for minimum E (eccentric anomaly); in 
other words, transfers were achieved to the desired orbit in the least amount 
of revolutions. The significance of this type of transfer trajectory is that 
it lends considerable insight into the minimum time, minimum fuel problem. 

The problem formulation involves deriving the perturbation equations 
of motion for fixed thrust engines of low, constant magnitude, and then opti- 
mizing them for minimum time and minimum E using Hamilton , Lagrange, and 
Mayer techniques in the calculus of variations. (Reference 4). After being 
optimized for minimum time and minimum E respectively, the differential 
equations are set up to solve for the optimum thrust angle necessary at all 
times to achieve the optimized transfers. These sets of equations are then solved 
by numerical integration,in this case on an IBM 360-30 computer system. An 
analytic solution is not available at present because of the inability to solve 
for Langrangian multipliers. This is the reason for the numerical approach. 

The basic assumptions follows: 

1) The thrust acceleration is much less than the local gravitational acceler- 
ation. (A<¢cg). 

2) Continuous, fixed-magnitude thrust. 

3) Relatively small total fuel consumption; therefore, for the sake of simpli- 
city, the thrust acceleration is assumed to be essentially constant over the 
transfer. That iS, constant thrust magnitude is taken aS synonymous with 
constant thrust acceleration. (A program modification to incorporate the 
changes in thrust acceleration due to fuel depletion is simple to implement). 


A description of the problem and its formulation follows. 





II. DESCRIPTION OF PROBLEM 

An initial elliptic orbit is assumed about a planet or any other 
body which has an inverse-square gravitational field. The ellipse has 
specifications (ao> eo) > where a. is the initial semi-major axis and 
ey is the initial eccentricity. The problem then is to transfer optimally 
to a desired final coaxial, coplanar ellipse with specifications (a,, ec)- 

See figure 1. This is accomplished by turning on the engine at the initial 
ellipse periapse and letting the optimally-directed fixed thrust program guide 
the vehicle along a spiral trajectory until the final ellipse is obtained. 

The engine is then shut off and the spacecraft is in its new orbit, the trans- 
fer having been conducted in minimum time and fuel expenditure or in the least 
amount of revolutions, depending upon which program is used. The numerical 
results show that the minimum time transfers always take more revolutions 

than the minimum E transfers, and the minimum E transfers always take more 
time than the minimum time transfers. This is as expected and lends consider- 
able support to the calculations. 

Each particular transfer is determined by the selection of the initial 
Lagrangian multiplers. Given an initial ellipse, different sets of inital 
Lagrangian multipliers will cause the vehicle to transfer to different ellipses. 
Since we have no analytic way of picking these multipliers at present, many 


cases are run to cover a wide range of transfers from given initial orbits. 


The formulation of both the minimum time and minimum E equations 


follows as the first phase of the problem solution. 
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FIGURE 1 
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tangential thrust forces, S and T respectively, shown in Figure 2, 


LAGRANGE FORMULATION IN VARIABLE E 
Orbital Perturbation Equations 


The general orbital perturbation equations in terms of radial 


are as follows: (Reference 5) 


Some useful relations follow: 


a z 
T > Me 
ao” « 
\ a 
oa 
Oo” 
m (VEHICLE) 
M 4 
FIGURE 2. 
Vehicle Force Vectors 
: daw. 2 ‘ =p) 
aes pe—, (Sesinf +T—) 
nVl-e 
J 2 
e = ce ssinf + T(cos E + cos f) 
dt na 


Using: 


1 ee 
“a ee ) 


cos E + cos f£ = 


(1) 


(2) 











p = a(l-e°) 
f = true anomaly = 0 
E = eccentric anomaly 


Equations (1) and (2) become; 


2 

n aC = Y= » = ,k = GCM 
3 

= k = 
. ue OC Gpecilody 
S = FNORMAL/ m = F sina / 
T= F PANG. /m F cosa/a@ 

For an ellipse; 
ee lees) 
l+e cos f 


Through substitution and simplification equation (3) yields: 


1 Fs [ae . | 
= = ~ Fag S : | sin 8 + sina + (1 + e cos 8) cos 


(3) 


(4) 


(>) 


(6) 


@) 


(8) 


(9) 


(10) 





Similarly, equation (4) yields; 


2 
che Ie NE Zac Osmortie (04 COS) 10) (COsia 
e = ce ae Nae i Q@ sina + We coen0 (ial) 


Next, it is desired to write equations (10) and (11) in terms 
of the eccentric anomaly E. The following relations between E and @ 


are used (Reference 6): 


sin og = sin E (12) 


l-e 
1+ ecos 0 = Ie ese F (13) 


ise sin E (14) 


- sin @ = l-e cos E 
Be 6 ie cos E - e 
l-e cos E (15) 
P l- 4 
Cs aGs2 ae a(l-e cos E) (16) 
lt+te cos 0 i 
l-e 
l-e cos E 
1+ecos 9 = ea (17) 
1 - ecos E—E 
: 1g 3 
Sin E ay sin 0 (18) 
ae A 
d@ = ~—————._ dE (19) 





= Ole 


B. Comparison with Edelbaum's Equations 
Substituting equations (12-19), equations (10) and (11) can be 


written in terms of E: 


M [om - 


(l-e cos E)m 


2F c Sina sinE + cose? (20) 


: Z Ny 
e = ae [ae ) sin E sina + tee (2 cos E - e poce E-e) cosa | 
(21) 


Equations (20) and (21) agree with Edelbaum's equations for 


low thrust perturbation (Reference 1). 


3 
dt = dE i (l-e cos E) (22) 


equations (20) and (21) become: 


Since 


3 
da _ 2fAa 2 
aaa - c sin E sina+WVl-e cos | (23) 


Zi 
se 42. fare’ sin E sina + ae (2 cos E-e cos E-e) con (24) 


Since we are optimizing the thrust angle for minimum flight time 
the following integral is minimized: 
E 


t. = '{ /. (l-e cos E) dE (25) 
E 
O 





SIOe 


This is of the form (Reference 4, 7) 


E 
f 
fe f* (X)> X5, EB) dE (26) 
E 
O 
where: 
foo 5 (l- os E) 27 
we B =Wa ly. exc (27) 
X, = a5 X, =e (28) 


C. Employing the Hamiltonian 
To find the minimum extremal, we write the Hamiltonian for this system of 


equations, 


2 
H=-f +) dd, (CE) £, ; (29) 
O : i i 
i=l 
where 
dx 
1 da 
Sf ee a ee 30 
see dE dE e”) 
< - 2 Lge 
2 dE dE ; (31) 
The independent variable is the eccentric anomaly E. Substituting 
into H: 


d 
H = - ye"), (l-e cos E) + d << = ae (32) 





-ll- 


: : da de 
Substituting for dE and dE: 





| ay Dae 
He se) a / (l-e cos E) + r, L c Sin E sina + 


5) 2 
2 
ie. GOs - + de 7 jo ) sin E sina + 


Se (2 cos E-e font EB = 7e) cos (33) 


D. Euler-Lagrange Equations for Adjoint Variables 





The following are the Euler-Lagrange equations for the adjoint variables )_ 
t 


dx 
LL sO a 
dE = OX, (i V2) (34) 


Expanding equation (34): 

















dx of gf 
a i) 2 1 aB 
dE rr Sab a ag “ Qa ) + ( Nai = r ofa) a (39) 
da oe Ae 0B 
e l 2 1 
aes —— — f 
dE CxS ge gene ae) +, Cn i ove? ee (36) 


E. Transversality Condition 


Next, the transversality condition, 


f 
j HdE + A, x, = 0 (37) 
1 L 





oe 


which is to be evaluated using the boundary conditions: 


t = 0: a=a, e=e, E =0 (38) 


tJ 
ll 


OPEN (39) 


The choice EG = 0 is made for convenience; this is the periapse point, the 
usual choice for the initial retro-impulse into a bound planetocentric elliptic 
orbit from an incoming hyperbolic trajectory. 


Combining equations (37), (38), and (39), 


Aw) + Af, Siena ere 0 (at t=t -) (40) 
F. Control Equation for Optimal Thrust Angle 
Equation (40) is a necessary condition for a transfer to be optimized for 


minimum time (and fuel). For the optimal thrust angle a , 


OH _ : eau meee 41 
eee. oe tT fe ae (41) 
which becomes 
A 2A 2 
0 = a [e sin E cos@ - l-e sina] + 


Zz 
A Aa 


. ? 
— [(1-e) Sin E coS@ - yi-e (2"€05 E = encos Siteoee) sina] (42) 





Se 


Divide through by cosa@ and solve for tana 


2a Xv Pe Sin E + ae (1 a4) sin E 


tan a ee eerie rere (43) 
2 Bos 
are Vl-e + Ae yi-e" (2 cos E- e oe E =e) 
To simplify the equation for computer programming purposes let: 
X = esinE (44) 
2 : 
Zeeazee -e)) Sin E (45) 
s 2 
ee l-e (46) 
a \ Z 2 
v= l-e (2 cos E - e cos E - e) (47) 
Neda KX Z (48) 
a e 
Dios 2a X% Y+ A.W (49) 
a e 
Substituting into equation (41): 
2a A - xe A Z N 
ere a Yt AW OD — 


This is the control equation for the optimal thrust angle a. The step variable 
of integration is E. Time of transfer (t) is determined by integration of 


equation (22) which can be written; 


co ee 


wee SS 
— 


i 
dE 8B (515 





Zien 


Thus, we now have sufficient equations for numerical solutions for 


different transfers, depending on the values assigned to Nis and x 
e ° 
O O 


The problem is next solved using the Mayer formulation to develop an 


alternate but equivalent set of equations. 


—— 
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IV. MAYER FORMULATION IN VARIABLE t 


The Mayer formulation of the problem proceeds as follows: 


dx. 

Cm ice Ge G2) 
Where: 

X, = a, X, =e, X, = E, Bee (53) 


The independent variable now is t rather than E as in the 


Lagrange formulation; here we minimize the value of t., which will be 


1 
called G. Thus; 
= 4 

G=t, (35) 
A. Hamiltonian and Euler-Lagrange Equations 
The Hamiltonian is different and does not contain f. 

3 
H= »), AGO) (55) 


i=l 


The Euler-Lagrange equations for the adjoint variables dh; are 


dA, 

















iL 0H 
ee 36 
dt OX, ca 
i 
Thus; 
dy of af 
a@ == 1 1 2 29 af (57) 
dt wa Qa aS pa) 72 Caf tAcl2 * Ae) Va 





=ihes 














d x of Aue 
2. ale 1 Z 1 of 
dt B Se de ae oe ee po a x 12 de (2g) 
da af af 
E 1 1 2 1 of 
= ll er ane f ae 
dt a ‘*a a * *e gr ?* A (12 Ae ee (59) 
B. Control Equation for Optimal Thrust Angle 
For the optimal thrust angle: 
OH = 0 = xr ory + oO) 
Oa ada € da (60) 


Equations (57), (58), and (60) are identical to equations (35), (36), and 
(41), using dE = 7 dt. Thus, the equation for tana is the same 
as equation (50): 
tan Benga ke oe (61) 
2 2a rx Y+ AW - D 
a e 
Also, since t is not explicit in equation (52), the Hamiltonian 


1s a constant: 


q = ——___ = CONSTANT (62) 


C. Transversality Condition 


The transversality condition is now 


f 
| - Hdt + A ox] + dG = 0 (63a) 
O 


There fore f 


. Hdts+oN dats) de +a, 0 | + dt = 0 (63b) 
a e E f 





=16a- 
Applying the boundary conditions from equations (38) and (39), 


dt; da; de; dE da-; de, = 0 (63c) 
: + + = 
Hdt . AE¢ dE dt . 0 (63d) 
Hence, since Ee is unspecified, 
and since te is also unspecified, 
- Hdt , ap dt . = 0 (65) 
which yields H = 1 (66) 
Combining equations (62), (64) and (66) 
7c! a ae mane | 
Delos Sea (all t) (67) 
B 
and: 
= 4,4) * 46% = fey 0, (t = t,) (68) 


Using equations (67) and (68) the Euler-Lagrange system of three 


equations (57), (58) and (59) can be reduced to the following two equations: 


Ot a) o 
The solution of equations 
ac z (71) 
c- “! (72) 
ee - 7 (73) 


together with equations (68), (69) and (70) is obtained numerically by choosing 
initial values as follows: 
tC = = 
BAe iy HE Ee AO) oy oc) (74) 
This requires guessing of values for A (0) and A 09) since there is 


at present no analytic way of choosing them. They must be chosen so that the 


humerical integration yields the required set of final values as follows: 


a=a., ee, i, = 0 (75) 
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V. THE MINIMUM E PROBLEM 
As mentioned earlier, the minimum’ E_ problem sheds some light 


on the minimum t problem. We now will minimize the final FE (number of 


turns goi 210) a = 
going m aos eo» ES O to ae > e,)- 


A. Hamiltonian and Euler-Lagrange Equations 


Using equations (30) and (31) minimize G = Eee The Hamiltonian 


is 
n= r a‘l + dA afo (76) 


The Euler-Lagrange equations are 

















da af of 
pee OH 1 2 
dE i da ~ ( Aa da : eran ) Ce 
ove =< _OH = = ( r an + r i: ) 
dE de a oe e oe (78) 


These are different from the Euler-Lagrange equations for minimum 


t optimization. 


B. Control Equation for Optimal Thrust Angle 


For the optimum thrust angle: 


OH af, df, (79) 











Zales 


C. Transversality Condition 


Now, writing the transversality equation, 


i: 
va = X on + dG 
ieee 


O 


where 


and expanding equation (80a) 


O 
-Hdf - HdE 
O £ 


O 
» ad + dE, = 0 


There fore 


- HdE. + dE, = 0 


which yields 


O O 
i \ da + rad + add + 


(80a) 


(80b) 


(80c) 


(80d) 


(80e) 
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VI. TRAJECTORY ANALYSIS 
A. Comparison of Minimum t and Minimum E Trajectories 
Equation (80e) is significant; it means that any pair of initial 
A's Qa, » Xe) will automatically generate minimum E_ trajectories 
terminating at any Ae» Of ; and E » This is so because the condition 
of equation (80e) can always be met simply by scaling the A's . Note that 
equations (77) and (78) are both linear in the X's and that the optimal 
thust angle control equation (79) relies only on the ratio of A. and 
not on their separate magnitudes. Automatic satisfaction of equation (80e) 
makes it easy to calculate minimum E_ trajectories and to construct a 
graphical representation of minimum E transfers on an e/a _ plane. 
Consequently, many transfer trajectories are generated which optimize 
E , and also some which optimize t . These optimal orbital transfers can 
be seen in Figure 3, which is a graph of e vs a/a . . It shows both minimum t 
and minimum E trajectories which are generated from the same initial orbit 
around Venus. These particular transfers concern a problem incorporating 
solar-electric orbiters around the planet Venus. (The thrust acceleration 
used in the calculations is oon ae The initial gravitational acceler- 
ation is about 0.8 m/s”; thus the basic assumption that A¢€g_ is satisfied. 
In an actual case the value of thrust acceleration would be about 
25: of 5x107> Vice However, the value of 5x10” is used for numerical 
calculations in order to speed up the trajectory and computer times, and thus 


reduce computer costs. Thus all calculated times must be scaled up by a 


factor of 25 to achieve a more realistic mission time. 





290 


Figure 4 shows minimum E trajectories only, with lines of 
constant t and E superimposed. Examination of the minimum t_ and 
minimum E trajectories readily discloses that orbital transfers from one 
point on the graph to another are accomplished in the least time by following 
the minimum t_ path and in the least total revolutions by following the 
minimum E path. Figure 5 shows the minimum t trajectories superimposed 
on the constant t and E lines of Figure 4. The minimum E_ trajectories 
are omitted for clarity's sake. 

The actual calculated orbital elements oscillate over each turn, 
as illustrated by curve 1A in Figure 5. For simplicity, the remaining curves 
are constructed by drawing a smooth curve through the values of the orbital 
elements of each turn at E = 2nz7 (n = 0, 1, 2, ...). Curve 1A shows the exact 
path of a minimum t trajectory. It is obviously an oscillatory type motion, 
which is much more pronounced here than it would normally be due to the 
increased thrust acceleration (25 times greater than normal), With the actual 


4 


thrust of about 2x10 ee the oscillations would adhere much more closely 

to the average path depicted by the heavy dotted lines - (curve 1B) in 

Figure 5. The minimum t trajectories, which are also minimum fuel trajectories 

here, are very similar to the minimum fuel trajectories calculated with variable 
Ke 

thrust by Edelbaum. 

B. Rotation of Line of Apsides 


When considering optimal orbital transfers between coaxial orbits, one 


wishes to maintain the same line of apsides. Since this program is designed 


OS ee 8 ee Ee Oe 
* 
See Section VIII. 
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to change only a and e =, no change in the line of apsides is expected. 

A close watch is kept of the change in apsidal line angle throughout all 
orbital transfers and it is verified through calculated results that such 
changes are small and relatively insignificant. In some instances a change of 
up to one degree in the line of apsides is detected over a transfer involving 
up to about 50 revolutions, i.e., one degree change in 18,000 degrees of 
revolution. 

Several points along curves 1B and 2 in Figure 5 indicated along 
with the flight time and eccentric anomaly elapsed at each position. 
Discussion of these points follows. 

Point J. Following a minimum time trajectory, the vehicle arrives 
at point J after 3.03 days' time and 37.69 radians of revolution. This 
point is clearly well past the three-day minimum E_ line and by inter- 
polation it is easy to see that a minimum £ trajectory would take about 
3.2 days to transfer to the same orbit. Now, while the minimum time 
trajectory has gotten to the same point in less time than the minimum E 
trajectory, it has gone through more degrees of revolution to accomplish 
this. Thus the minimum time trajectory has gone through 37.69 radians of 
revolution; while, by interpolation, it can be seen that the minimum E 
trajectory will reach the same point with only about 36.8 radians of 
revolution. In summary, a minimum time transfer takes 3.03 days in 37.69 
radians, while the minimum E transfer takes about 3.2 days in 36.8 radians. 


Point L. Here the minimum ct trajectory has reached che point L 
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in 4.57 days while it takes the minimum E trajectory about 4.8 days. How- 

ever, it takes the minimum t_ trajectory 50.26 radians of revolution as 

opposed to only about 47 radians of revolution for the minimum E _ trajectory. 
Point P. This is a very good point for comparison. The minimum t 

flight time is 4.9 days. This is clearly less than the minimum E time, 

which is well past the five day line. The eccentric anomaly elapsed is 

46.9 radians which clearly exceeds the value for the minimum E _ trajectory 


since point P is distinctly below the minimum E_ contour line where 


E = 45.23 radians. 
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VII. DISCUSSION OF OSCILLATORY MOTION 


A. Effect of Thrust Acceleration and Lagrangian Multipliers on Oscillatory Motion 

Through numerous computer runs, it was verified that different 
thrust accelerations produce essentially the same mean transfer trajectories, 
providing the Lagrangian multipliers are scaled according to the ratio of 
the accelerations (observe Figure 6). A very low thrust system produces a 
trajectory with minimal oscillatory motion (curve Ay) + As the thrust 
acceleration is increased, and the Lagrangian multipliers scaled accordingly, 
the oscillatory motion of the transfer trajectory becomes much more 
pronounced (curve A,)> while the mean trajectory remains the same as with the 
lower thrust acceleration. 

It has also been verified (Figure 7) that minute changes in the initial 
Lagrangian multipliers result in uniquely different transfer trajectories 
which do not overlap, as long as the same thrust acceleration is used. The 
transfer trajectories can overlap if both the Lagrangian multipliers and the 
thrust acceleration are changed (Figure 8). 

Figure 9 shows orthogonality between minimum E trajectories and 
constant E path which must exist for optimal solutions. This is a good 
check on the problem solution. In this case, each point was plotted to show 
the detailed angular orientation of the two curves at the point of inter- 
section. The smoothed-in path of trajectory P in Figure 4 does not appear 
to show orthogonality, but close observation, i.e., blow-up in Figure 9, proves 
it to be so. 

B. Appearance of Negative Eccentricities 


An interesting and perplexing problem was encountered while examining 
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the oscillatory points along trajectories at very small eccentricities. 

The problem is one which involves explaining the presence of a few negative 
eccentricity values which appear during the calculation of transfers of low 
eccentricity. Figures 10 and 11 illustrate this point. In the blow-up of 
area I in Figure 11, points A, B, C and D all lie below the axis and are 
thus negative eccentricities (of very low magnitude) and are printed out as 
such by the computer. 

These negative points appear to be in line with the normal cycle 
of the oscillatory motion. Another point to note is that none of the 
smoothed trajectories ever go below zero eccentricity; that is to say, at 
the end of every cycle (i.e., for every Dae of eccentric anomaly) e 
is always greater than or equal to zero. However, to date there is no clear 
meaning accorded to these negative eccentricities. Further analysis is 


necessary to determine the real significance of this factor. 


C. Oscillatory Motion and the £ Function 

Also connected with the oscillatory motion which accompanies all 
spiral trajectories is the &£ function as defined in equation (40). Through 
many minimum t program runs, it was found that the € function never Stays 
at zero but always oscillates through zero once every revolution, usually at 
every 2nz of eccentric anomaly. The smoothed €& function is zero, but the 
minimum t formulation states that the € function must equal zero at the 
temninal point in order to represent an optimal transfer. Therefore, the 
transfers are exactly optimal once every revolution when the € function is 
zero and are near optimal in between. The minimum E_ transfers are exactly 
optimal at all points along the trajectory since there is no &€ function 


requirement in that formulation. 
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FIGURE 11 
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VIII. COMPUTATIONAL COMPARISON WITH EDELBAUM'S VARIABLE THRUST CALCULATIONS 


A. Examining Thrust Angle Schemes for Dif ferent Trajectories 


An optimal variable thrust system can be utilized to measure the 
limit of best possible performance with an optimal fixed-thrust system. 
This is so because of the extra degree of freedom added to the control system 
wnen the variable thrust is incorporated, making it a more efficient system 
than the fixed-thrust one. Thus, any optimal performance by the fixed-thrust 
system can, at best, be only equal to the optimal variable thrust system 
and, more often than not, it will be less than the variable thrust system's 


performance. 


The following equations were derived by Edelbaum (References 1, 2, 3) 


for optimal variable thrust trajectories using an averaging technique and are 


2) 


used here for comparison with the fixed-thrust performance: 


. —_— a iP 
, | - -1 
V eet. = A Ct 2 = = | tae 2 2 “= cos| 7 { sin ! e - sin 


c if f£ a 
O 


4) 


== Z 2 
A Z, ea 
(av .)* sat Tt = (2) see fae)” 


a’ = the averaged squared variable thrust acceleration over a 


- turn = a constant 
t 
i ——— 
2 2 
anos AZ — AL 
J= 9 dt 9 C, 


(81) 


(82) 


(83) 


(84) 
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tS (85a) 
, A ap eee 

P = POWER = 5 | | v, (85b) 


Equation (81) applies to large changes in (a, e); i.e., more 
than one revolution. Equation (82) is a special case of (81) and applies 
only to changes within one turn. 

To compare the two propulsion systems, Ne and AV. for fixed 


thrust are both defined as: 


Vv = AV. = Act (AC) (86) 


As can be seen in Table 1, the fixed thrust _ and Es are 


always larger than the corresponding variable thrust - or av. where: 


= _(4V)_ fixed __ 
c= AV) fixed (87) 


1~ (AV) variable 
Note that C, is always greater than 1. This indicates a higher cost for 
the fixed thrust system as expected. 

As might be anticipated in cases involving near-circular orbits 
the thrust magnitude of the variable thruster is almost constant, yielding 
the same case as the fixed thrust system. Here the value for C, is very 
close to one. When transfers are conducted with more elliptical orbits, 
the variable thrust system thrust magnitude scheme no longer resembles 


the fixed thrust case, and the cost is then less, making C, take on larger 


values. 
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B. Maximum Aa and Ae Steering Programs 


It is interesting and instructive to see how the thrust angle varies 
throughout each transfer. Though only the angle varies (and not the magnitude 
of the thrust) in fixed thrust transfers the angle program agrees closely to 
the one calculated by Edelbaum (See Figure 12). The thrust at points A through 
E coincide in both direction and magnitude. 

Examining single turn characteristics, if e is changing more 
rapidly than a, then the second term on the right-hand side of equation (82) 
will dominate the first term. Similarly, if a is changing more rapidly than 
e , then the first term dominates the second term. Figure 12 shows the 
thrust schemes for maximum Aa/a and maximum Ae programs. The last two 
columns in Table 1 show that when the predominant change is in a ,_ then the 
thrust scheme is an '"a'' thrust scheme and vice versa. 

The fact that Table 1 indicates close but slightly inferior 
performance for the fixed thrust when compared with the variable thrust lends 


considerable support to the validity of both sets of calculations. 
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IX. DETAILED ANALYSIS OF THE THRUST ANGLE SCHEME FOR A MINIMUM E 
TRAJECTORY 


The minimum E trajectory L from Figure 4 was selected for 
examination of its thrust angle scheme. The optimal, computed thrust angle 
was recorded and plotted versus eccentric anomaly (see Figure 13). The 
trust angle positions were recorded every 90 degrees of eccentric anomaly 
for every third revolution. 

Figure 4 shows that for the first 2.6 days the thrust program 
changed primarily the semi-major axis. In this time period a doubled from 
ZO Tee ow meters while the eccentricity increased very slowly from 0 to 
about 0.03. This part of the transfer is depicted in Figure 14 where the 
original orbit is 0, , and after 2.6 days the trajectory has reached 


1 


orbit 0, . Thus the program starts out by increasing the semi-major axis 
with little change in eccentricity; it uses essentially a maximum Aa_ thrust 
angle scheme with mostly tangential thrust. This can be seen by observing 
the thrust angle scheme for revolutions one and four in Figure 13. 

By observing trajectory L in Figure 4 one can see that between 
t = 2.6 and t = about 4.6 (days) the thrust program changes to a maximum 
Ae trajectory, while remaining fairly static in its semi-major axis. 
Figure 15 depicts this leg of the transfer from orbit 0. to orbit 0, 
These orbits were chosen at appropriate points to illustrate the computed 
trajectory. A drastic change in the thrust angle scheme between revolutions 
four and seven can be seen in Figure 13. 


Figure 16 shows the transfer between t = 4.6 and t = about 11 days. 
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The semi-major axis a decreases steadily while e increases. The engine 
is cut off when the trajectory reacnes orbit 0, . At this p@int e = 0799¢ 
and a= 1.24x10 meters. It is interesting to note that between t = 4.6 
and 11 days, as the progran begins to change a more rapidly, the thrust 


angle pattern (between revolutions 10 and 27) begins to revert back again 


toward a tangential type thrust scheme (Figure 13). 
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STARTING FROM ORBIT O, THE SPIRALING PROCESS HAS 
REACHED ORBIT O, ATt=2.6 DAYS. 
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AT t = 4.6 DAYS THE ORBIT SPIRALING PROCESS HAS 
CARRIED THE VEHICLE TO Oz. 


SECOND LEG OF MINIMUM E TRAYECTORY L 
SHOWING MAXIMUM Ae STEERING SCHEME 


FIGURE 15 
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AT t=Il DAYS, THE VEHICLE HAS BEEN 
TRANSFERREO TO ORBIT O, AT THIS POINT THE ENGINE IS 
CUT OFF AND THE VEHICLE REMAINS IN THE SPECIFIED ORBIT O,g- 


FINAL STAGE OF MINIMUM E TRAJECTORY L 
SHOWING ATTAINMENT OF DESIRED FINAL ORBIT (O04) 


FIGURE1 6 
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X. MARS ORBITAL TRANSFER 
A. Description of Problem 

In order to gain a practical comparison with single impulse 
chemical retro transfers, realistic values were assigned to solar electric and 
chemical rocket engines to be used in an orbital transfer about the planet 
Mars. An initial elliptical orbit was chosen with large eccentricity and 
large semi-major axis (when compared to the radius of Mars). The problem 
involved transferring from this initial orbit into a much smaller ieee 
orbit with radius of 1.1 that of the radius of Mars. The periapse of the 
initial elliptical orbit is equal to the radius of the desired final circular 
orbit. 

The low-thrust spiral transfer is accomplished by turning on the 
solar electric engine at the initial ellipse periapse, causing the space 
vehicle to slowly spiral into the desired final circular orbit using the 
optimized thrust angle program for fixed thrust. The retro transfer is 
accomplished by using a single chemical rocket retro at the ellipse periapse 
causing the space vehicle to transfer directly into the circular orbit. The 
amount of fuel required for each transfer is then compared. The purpose of the 
spiral transfer is to deliver a larger payload to the final orbit. One must 
then decide if the increase in payload delivered is worth the extra time required 
to make the transfer. Figure 1/7 shows the two methods of transfer being 
compared, the upper drawing being the spiral transfer and the lower one depicts 


the single chemical retro impulse transfer at the ellipse periapse. 
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B. Initial Data Used 


The values used in this comparison were as follows (Reference 8): 


Solar Electric Engine 
6.597560. 


F = 
ae 
vi = 39,464.9 m/s 
Chemical Rocket Engine 
MG = 2890.49 m/s 
Both Engines 
Mou = initial mass in ellipse = 2049.5093 kg 


radius of Mars: 3400 km 





ky “3 4.293x10!° m/s” 
ars 
a 7 
anil 9.734764x10 om 
ep 
e=1- = 0.95773909 

a 
ell 


C. Chemical Rocket Transfer 


The Mathematical Development follows: 
Ye, = velocity at ellipse periapse 


E ok -< (88) 
Dp 


4.6389x10> m/s 


I 


< 
| 


= velocity at final circular orbit 


\: = 3,390 m/s 


AV =v “ a = 1,248.9 m/s (89) 


I 
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we mass in final circular orbit after 
chemical retro 


me, /m_ a oe Aviv: = 97 1248.9/2890.49 
-.4 
= e ae 365.2 
So Bee = (.652) (2049.5093) = 1,340 kgs 


Ame, = 2049 .5093 
-1340.0000 
709.5093 kgs 


Thus, the direct chemical rocket insertion requires 709.5093 kg 


of fuel with an initial total vehicle mass of 2049.5093 kg. This means 


that a total mass of 1,340 kg is delivered into final orbit by this method. 


D. Optimized Fixed-Thrust Acceleration Spiral Transfer 


Next, the same transfer is accomplished using the solar electric 
spiral technique. Assuming the same initial ellipse and total vehicle 


mass, the initial thrust acceleration is determined: 


A= F910) 


ell 


4 ie 





m/s 


Although assumed a constant, A will actually increase as the 


vehicle mass decreases due to loss of fuel. A mean value of about 3x107¢ m/s 


was used in the programs, which implies a total anticipated change of about 


6 percent in thrust acceleration. The time to spiral optimally to Ce 


a 


(90) 


2 


f 


was computed to be 128.685 days. To compute the fuel expended while spiraling 


into the final circular orbit (Reference 9): 


F =|nr 

| 1h | v, 
lal ~ F. _ -59756046 een ee 
a v, 39,464.9 m/s 


is setae oar 


=e 
Hl 


(91) 





SHo- 


To find the fuel expended while spiraling, 4m, » multiple | a | 


by the totaltime t elapsed during the spiraling process: 


am. = J] e (92) 
Now, 
C =923-635 = {11 7terlon sec 
_ => 7 
Am. = (1.515x10.~ ke/s) (1sll716x10. ‘s) 


Ke 
Shoe ke 
E. Results and Comparison of Data 


Thus, the solar electric spiraling process has accomplished the 


same transfer while using only 169 kg of fuel. This yields 


Initial mass 2049 .5093 
Fuel - 169.0000 
ce 1880.5093 kg 


1880.5093 kg of mass is delivered into final circular orbit by the 
solar electric method. Since 709.5093 kg of fuel were required by the chemical 
retro transfer, and only 169 kg were expended during the solar electric spiral, 
the latter method yields a savings of 540.5093 kg or about 1150 pounds. 

These results mean that the desired final circular orbit can be 
obtained almost instantaneously by chemical insertion, but with a cost of 
709.5093 kg of fuel; or, the same orbital transfer can be accomplished over a 
period of 128.685 days with a cost of only 169 kg of fuel. Thus, if one chose 
the slower spiral technique, the amount of payload delivered to final orbit 


could be increased by about 1150 pounds. The advantages and disadvantages of 


* ° 
Note that for constant thrust magnitude, the total increase in thrust 
acceleration due to fuel depletion = 9 percent. 
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each method must be weighed along with the mission objectives to arrive at 

a decision over which method should be used. It would seem that the spiral 
technique as described would be far better suited to an unmanned observational 
missions than the chemical method (Reference 10). It would seem impractical 
to use the slow spiral technique on any manned flights, regardless of the fuel 
saved, but this does not preclude the use of solar electric engines enroute 


to various destinations as a thruster boost. 
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XI. CONCLUSIONS AND RECOMMENDATIONS 

The numerical solution of the optimum fixed-thrust acceleration 
elliptic spiral technique has been used to generate many different planeto- 
centric transfer trajectories. Figure 2 represents an infinite number of 
minimum time and minimum eccentric anomaly transfers. It should be noted that 
any transfer generated can also be reversed by appropriate manipulation of 
the Lagrangian multipliers. Optimized transfers can be conducted from any 
point along a generated trajectory to another point on that trajectory. This 
is done by starting the program with the Lagrangian multipliers already 
printed out at the desired point along the trajectory. 

At present there is no analytic method for choosing the Lagrangian 
multipliers for particular transfers. By generating many transfer trajectories, 
one can accumulate many cases, but the initial process involved in choosing 
the Lagrangian multipliers remains guesswork and experimentation. Possibly 
by examining more closely the transfers produced by the numerical solutions 
used herein one might draw some clues for the analytic solution. The use of 
reiteration techniques to converge initial guesses to final correct values 
was not explored in this thesis, but is a well-established art (Reference 11). 

As mentioned in an earlier section, all transfers were designed to 
be coplanar, and co-axial. As expected, virtually no rotation of the line of 
apsides was observed. A problem noted and recommended for further study was 
the occasional production of very small negative eccentricities along certain 
transfer trajectories. It only occurs in transfers involving circular or 


near-circular orbits and then only in transfers displaying maximum 4a thrust 
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angle schemes. 

The hypothetical Mars orbital transfer conducted under realistic 
cond it ions yields some interesting facts for consideration. The optimized 
spiral technique proves to be able to deliver up to 1150 pounds more into the 
desired final orbit than the single impulse chemical retro method, based upon 
a decrease in the amount of fuel required to accomplish the transfer. However, 
the spiral See requires 128 days to make the transfer, while the chemical 
method works almost instantaneously. 

Another possibility would be to assume that no more payload is 
required on this particular flight, in which case the extra 1150 pounds 
could be trimmed from the initial mission weight. This would greatly reduce 
the fuel required to launch the vehicle from earth orbit; or, one could use 
the same amount of booster fuel in which case the vehicle would reach the 
target planet much faster than before and would erereny decrease the overall 
flight time. One must then look at all of the factors involved in particular 
missions to determine which method or combination of methods might be most 
useful. 

While this thesis has addressed itself only to the use of optimized 
fixed-thrust orbital transfers, it is quite conceivable, and possibly even 
more advantageous to consider using the electrostatic engine in heliocentric 
Space enroute to the target planet as well as after planetocentric space has 
been entered. Any judgment on how to interface a heliocentric interplanetary 
trajectory with planetary orbit insertion will have to be based upon further 


study which must involve optimizing the entire flight from earth orbit to the 
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desired final orbit about the target body. Figure 18 illustrates this type 
of problem. This involves the incoming parabolic trajectory, chemical retro- 
insertion into a very large eccentric ellipse, followed by the optimized 
spiraling transfer technique until the desired final orbit is obtained; in 
this case the final orbit is a small circular one. 

It seems obvious that the slow spiral method of transferral would 
not be desirable for manned flights because of their long duration, but 
there are some instances where such an alternative might be preferred. For 
unmanned exploratory flights, where time may not be the most significant 
factor, this method offers an excellent chance to study the upper atmosphere 
of a given planet and to make planetary observations in a low, systematic 
manner. Starting in a large elliptical orbit, a spiral trajectory which 
steadily decreases in altitude until a low circular orbit is achieved might 
be desirable to record radiation, magnetic, and electrical fields at 
systematically decreasing altitudes, along with other upper level atmospheric 
conditions, planetary observations, etc. Orbital changes conducted by large 
impulse chemical retros make it impossible to record such long-time data as 
the transfers occur so quickly. There are many factors to consider when 
weighing the merits and drawbacks of optimized electrostatic propulsion 
(Reference 10). 

The surface has only been scratched in this field. Improvements 
lie in many directions. For example, Handelsman (Reference 12) has put forth 
a formulation for the introduction of a switching function into the problem 


which would optimize transfer trajectories for minimum fuel and fixed-time 
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restraints. This would allow part of the trajectory to consist of coasting 
arcs in which the engine would be shut off by the switching function in order 
to comply with the fixed~-time restraint. 

Based upon the work done to date the following recommendations for 
further study are made. 

1. Investigation of the occasional occurrence of negative 
eccentricities in certain transfer trajectories. 

2. Use of a thrust acceleration which, for fixed-thrust magnitude, 
increases as the vehicle mass decreases due to fuel expenditure. 

3. Optimization of the entire flight trajectory including inter- 
facing with the rendezvous and final orbit insertion in order to find the best 
place to use electrostatic thrusting in the flight plan. 

4. Study the solutions generated with a switching function for 
fixed-time: problems and compare with solutions of the type covered in this 
thesis. 

5. Addition of other orbital changes, such as rotation of apsidal 
line, non-planar transfers, etc. 

6. Investigation of an analytic solution for a method of selecting 
the correct initial adjoint variables required for particular transfers 
desired. This would eliminate the guesswork now required at the outset of 


every problem. If this is not successful, iteration techniques may be 


necessary. 





Pl. 


12. 


13. 
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APPENDIX A 


Description of Computer Program 

This is a description of the computer program used for the numerical 
solution of optimal (minimum time and minimum E) planetocentric spiral 
trajectories using fixed power, fixed thrust propulsion, where the thrust 
acceleration is much less than the local gravitational acceleration. The 
program is written in Fortran IV computer language and incorporates a 
straight forward Runge-Kutta type integration process. It was designed for 
use on the 360/30 computer system, but can easily be adapted for use on 
other systems. Double precision was used in all program calculations 
(Reference 13). 

In this appendix the important program input and output parameters 


are described, followed by the program itself and the program flow chart. 


Program Input Parameters 


In order to produce a trajectory for optimal orbital transfer, 
certain initial quantities must be fed into the computer. These can be 
divided basically into three sets of values; the initial orbital conditions, 
the two Lagrangian variables Omen: and the thruster conditions. 

The initial orbital conditions include: 

a (initial semi-major axis) 
e (initial eccentricity) 


k (mass constant) 


t = 0 (transfer time) 
E = 0 (eccentric anomaly) 
O 
w =O (angular position of line of apsides) 





ez 


The thruster conditions include: 
A (thrust acceleration) 

Ml (number of revolutions desired) 

M2 (number of steps per revolution) 


The initial Lagrangian variables are: 


0)= 
A 6 JE dg 


A 
ef BA, 


The accuracy is determined by M2 (the number of steps per revolution 
desired). In other words, if M2 = 50, then the program will compute and 
print out the spacecraft's coordinates 50 times during each revolution, giving 
the exact description of the elliptical orbit which would be obtained were 
the engine to be cut off at any one of those points. However, we have found 
through experimentation that a value of M2 = 10 is sufficiently accurate and 
generates trajectories more quickly. 

For every initial set of Lagrangian variables, a different transfer 
trajectory will be generated. The type of transfers generated are indicated 


? 


in Figure 4. 


Computer Output Parameters 


At each point in the transfer trajectory a number of instantaneous 
quantities are printed out. These include the total step number, which, if 
divided by the number of steps per revolution, yields the number of revolutions 
to that point (a total step number of 156 would mean 15.6 revolutions at 


ten steps per revolution); the eccentric anomaly in radians; the thrust angle 
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in degrees; the & function, which must go through zero for an optimal transfer 
to take place; the eccentricity; the change in eccentricity per radian of 
eccentric anomaly; the semi-major axis; the change in semi-major axis per 
radian of eccentric anomaly; the Lagrangian variable r, and its change per 
radian; the Lagrangian variable 5 and its change per radian; the time 
elapsed since thrust initiation; the number of seconds of time per radian of 
eccentric anaomaly; angle of the line of apsides; the change in the line of 


apsides per radian of revolution. 
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APPENDIX B 


Program Listing 


Cava MAIN 
DOUBLE PRECISION W,X,Y,2,h4,GMU,N,D,H-SIGM 
COMMON W.X,7Y.Z, A, GU PN, D,H,SIGN 


DOUBLE PRECISION H1(4),H2(4),H3(4) ,0(7) ,00 (7) ,DQ(7) ,DELQ(7),D1, 02, 
rr. PRT) 

NAMELIST /INPUT/ M1,M2,A,GMU,0,SIGN 

DATA QO / 7*0.D0 / 

SIGN = 1.D0 


pO 10 M=1, 4 
H2¢e) = 1.D0 
H1(H) = 0.5D0 
10 43(M) = 0.5D0 
H1(3) = 1.D0 
U3(4) = 1.D0 
ioe1) = 0.D0 
H2(1) = H2(1) / 6.D0 
H2(4) = H2(1) 
Men2y = 1212) 7 3. D0 
H2(3) = H2(2) 


20 READ (1, INPUT, END=200) 
IF( H1) 25, 200, 25 
PoeniTTe (3,40) Ni, “2, A, GHU » SIGN 
4UQ FORMAT( 1H1, I5, * REVOLUTIONS", I5," STEPS PER REV. A = ', 
jetep15.8, © MU = *,D1S.8 ,* SIGH =',p15.8 } 
WRITE (3,INPUT) 
T = 0.DO 
KMAX = IABS( M1 * M2 ) 
O(7) = Q(1) 
X = M2 
DT = 6.2831853071795865 / X 
KOUNT = 0 
TOO TT = T 
DO 101 M = 1, 6 
PR (+1) = O(M) 
Pom) = O(%) 
one eLO(M) = 0.DO0 
feos) = 0. DO 


QQ(7) = Q(7) 

PR(1) = T * 57.2957795 
PR(G) = PR(6) / 86400.D0 
PR(7) = 0.D0 


IF (0 (6) -NE.0.D0.OR.Q(7) -NE.0.D0) PR(7) = DATAN2 (Q (6) ,Q(7) ) 
PR(7) = PR(7) * 57.2957795 
WRITE(3,50) KOUNT, T, PR 
50 FORMAT( 1HO, I6,1P8D14.6 ) 
pe peniy ( TO, -DO ve) 
PR(1) = DATAN2(N,D) * 57.2957795 
PR(2) = DOQ(1) #0(3) + DQ(2)*Q(4) - DO(5S) 
WRITE(3,51) PR(1), PR(2), (DQ(M) -M=1,6) 
S51 FORMAT( 1H , 8X, 1P8D14.6 ) 
KOUNT = KOUNT + 1 


102 K = 1 + K 
=o t pr * Hs (4) 





ee 


104 


105 
106 


Oy 


120 


200 


B=2 


G@RELSOERIV( T, 0, DO ) 

D1 = HI(K) * DF 

D2 = H2(K) * DT 

DO 104 m = 1, 7 

DELQ(M) = DELQ(M) + D2 * DQ(M) 

O(M) = OO(M) + DO(M) * D1 

IF( K - 4) 102, 105, 102 

DO 106 M = 1, 7 

Cy = COMM) © DELO (MN) 

IF( KOUNT — KMAX ) 100, 107, 100 

GULL DERIV( T, QO, DO ) 

DO 120 % = 1, 6 

PR(M+1) = Q(M) 

PR(1) = T * 57.2957795 

PR(6) = PR(6) / 86400.D0 

PR(7) = 0.DO 

IF (0 (6) .NE.0.DO.OR.Q(7) .NE.0.DO) PR(7) = DATAN2(Q(6),Q(7) 
PR(7) = PR(7) * 57.2957795 

WRITE (3,50) KOUNT, T, PR 

PR(1) = DATAN2( N,D) * 57.2957795 

PR(2) = DQ(1)*Q(3) + DQ(2)*0(4) - DOQ(5) 
WRITE(3,51) PR(1), PR(2), (DQ(M) ,M=1,6) 
SONTO 20 

RETURN 

END 
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Cave DERIV 

IMPLICIT RFAL*8 ( A-H, O-Z ) 

SUBROUTINE DERIV ( E, Q, DQ) 

DOUBLE PRECISION W,X,Y,Z,A,GMU,N,D,H,SIGN 

COMMON Wek, ¥e2rh, GuULN.D,H, Spee 

DOUBLE PRECISION E,Q(4) ,DQ (4) ,DSIN, DCOS, DSQRT, DABS, SINE,COSE, 

imeesler, T1, T2, T3 

SORTAU = DSORT( Q(2) A Guu } 

SINS = DSIN( Ff ) 

Gor = DCOS( F } 

COS2ZE = COSE ® COSE 


Geel DO — O(1) * 0O(1) 

Y = DSQRT( DABS( Z) ) 

Z= 2% * SINE 

X = 0(1) * SINE 

Wee 1S. (2. D0%% COSE -— Q(1) “¥ *(COSZ2E+ 1.D0)} 
N = 2.00 * Q(2) * Q(4) * X # Q(3) * 2 

D = 2.D0 * Q9(2) * Q(4) * Y + Q(3) * W 

H = DSORT( N*®N + D¥*D ) * SIGN 

eee X*N + Y*¥D 


T2 = Z*N + W#D 
j= WMO (2)8/ GMU / i 


Honeys = 2.D0 & TI O(2) £ 0(2) «© 71 

DO(1) = 13 2 0(2) "4 T2 

DO(5) = SORTAU * (1.D0 - Q(1) * COSE) 

DO(8) = -2.D0 * T3 * (3.D0 # O(2) * O(4) * T1 + O(3) Bae) 


fet 16. 5D0° ADO (5) 

DO(5) = DO(S) * QO(2) 

DO) = =T3mO( 2)" (2. DO#0 (2) 40 (4) MUISSINE = O(1)* DAa es 
03) *(2.DO0FK PN + (0 (1) FW/ (YPY) + Y*(1.D0 + COSZE Ar 
2 - SORTAU * Q(2) * COSE 


DQ(6) = T3*Q (2) * (Z*D/ (Y*Y) *.(2.D0-O (1) * (COserO (ae 
V Y*N*(COSE - Q(1) ) ) 
T3 = 0.DO 


IF (0 (6) ~NE.0.D0.OR.0 (7) .NE.0.D0)T3 = DATAN2(0(6),9(7) ) 
T1 = DSIN (T3) 

T2 = DCOS(T3) 

DOQ(7) = T2*DQ(1) - T1 * DO (6) 

DQ (6) T2*DQ(6) - T1*DQ(1) 

RETURN 

END 








BEGIN DO LOOP 
1OM= 1, 4 













H2 (M)} = 1.p0 
Hl (M) = 0.5p0 








H3 (M) = 0.5D0 


END OF DO 
LOOP ? 






Hil §¢3)> = “poe 







H3 (4) = 1.D0 







H3 (1) = 0.D0 






H2 (1) = 
H2 (1)/6.D0 









H2 (4) 






= HZ. (1) 






HZ (02) >= 
H2 (2)/3.D0 






H2 (3) = H2 (2) 





READ FROM DEV 1 
VIA FORMAT 
INPUT 







END OF DATA? 
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Program Flow Chart 


EXIT 
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WRITE TO DEV 3 
VIA FORMAT 
40 


FROM THE LIST 








NOTE 13 


EMST = MEM 2 ee. 
GMU, SIGN 





WRITE TO DEV 3 
VIA FORMAT 
INPUT 





T = 0.DO 


KMAX @ 
LABS (M1L*M2) 


Q 7 2Q QM) 





100 1? 


NOTE 18 


BEGIN DO LOOP 
10iImM = 1, 6 


02.01 


19 


PR (M * 1) © Q (mM) 
QQ (MY = Q (») 


101 20 


DELQ (mM! = 0.D) | 


} 





01.20 





NO 





END OF DO LOOP? 


7 & 


02 
DELQ (7) = 0.D0 













QQ (7) = Q (7) 






PR (1) *= 
T*57.2957795 







PR (6) = 
PR (6)/86400.D0 






PR (7) = 0.D0 











Q(6) .NE. 
O0.DO .OR. 
Q(7) .NE. 
0.DO 


FALSE 













04 
PR (7) = 
DATAN2(Q(6), Q(7)) 


05 


FR (7) = 
PR (7)*57.2957795 


06 


WRITE TO DEV 7 
VIA FORMAT 50 







NOTE O7 


LIST = KOUNT, T, PR 


DER IV 


(T, Q, DQ) 
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PRC). = 
DATAN2 (N, D) 
*57.2957795 


PR (2) © 
DQ (1)*Q (3) + 
DQ (2)*Q (4) - DQ (5) 


11 





WRITE TO DEV 3 
VIA FORMAT 51 
FROM THE LIST 







LIST = PR (1) , 
PR (2), (DQ(M), M ® 
1,6) 





13 


KOUNT © KOUNT + 1 





















K = ] +k 


T° Tt + DIS) 


15 


DERIV r 
(T, Q, DQ) | 


103 16 
Dl = HI (K)*DT 


D2 = H2 (KY*DT 





NOTE 17 
BEGIN DO LOOP 
104 m = 1, 7 
18 
DELQ (mM) © 


DELQ (M) + D2*DQ (M) 


19 


Q (M) © QQ (Mm) + 
DQ (M)*D1 


20 
END OF DO LOOP? 
ES 
21 
-/+) 
O 


105 NOTE 22 


BEGIN DO LOOP 
LOG. er" UN 


03.01 
106 23 


Q (M) = QQ (mM) + 
DELQ (M) 


3.01 
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O2G2 
ol 
END OF DO nO 
LOOP? 
. (4) 
02 
2/4) 
. @ 
100 
107 03 
DERIV 09 
(T,Q,DQ) 





PR(7) = 
DATAN2 (Q(6) Q(7)) 
NOTE 04 10 
BEGIN DO LOOP PR(7) = 
120 mM = 1, 6 PR(7)*57.2957795 
l 


05 1 
PR(M + 1) = Q(M) WRITE TO DEV 
06 
END OF DO NOTE 12 


LOOP? 





PR(1) «= 


DATAN2 (N,D) 
#57,2957795 


PR(2) 
DQ(1)*Q(3) + 
DQ(2)*Q(4) = 

DQ(5) 


PR(1) = 
T*57.2957795 


PR(6) = 
PR (6)/86400.D0 





PR (7) = 0.D0 





WRITE TU DEV 
5 
VIA FORMAT 
$1 
FROM THE LIST 







NOTE 15 


LIST = PR(1L) 
PR(Z), (DQ(M),M = 
L5) 





SUBROUTINE DERIV (E,Q,DQ) 


DERIV 


02.08* 


SQRTAU = 
DSQRT (Q(2)/GMU) 


SINE = DSIN (E) 
COSE = DCOS (E) 
COS2E = COSE*COSE 


z= 1.D0 - 
Q(1) *(1) 


Ye 
DSQRT (DABS (Z) 
Z © Z*SINE 


X = Q(1) *SINE 


We 
Y*(2.DO*COSE - 
Q(1)*(COS2E + 

1.D0)) 


N e 
2.DO*Q(2)*Q(4) 
*X + Q(3)*Z 


De. 
2.DO*Q(2)*Q(4) 
*Y + Q(3) *W 


H * DSQRT(N*N + 
D*D) *SIGN 


Tl = X*N + Y*D 
T2 = Z2*N + W*D 
T3 = A*Q(2)/GMU/H 


DQ (2) « 
- 2. DOATS*#Q (2) #Q(2) 
*TL 


DQ(1) = 
T3*(2)*T2 


DQ (5) 
SQRTAU*(1.DO - 
Q(1) *COSE) 


pq(4) = - 

2.D0*T3(3.DO*Q 

(2)*Q(4)*TL + 
Q(3)*T2) + 
1.5b0*DQ(5) 


DQ(5) © 
DQ (51*Q(2) 








FALSE 


0? 





DQ(3) = - 
T3*Q(2) 
*(2.DO*Q(2)*Q(4) 
*(N*®SINE - 
Q(1)*D/Y - 
Q(3)*(2.DOKXaN + 
(Q(1) *w/(Y*Y) + 
Y*(1.D0 + 
COS2E))*D)) - 
SQRTAU*Q (2) *COSE 













DQ(6) = 
T3*Q(2) 
*(Z*D/ (Y*Y) 
*(2.D0 - 
Q(1)*(COSE + 
Q(1l))) - 
Y*N*(COSE - 
Q(1))) 


T3 = 0.D0 





13 
DANTAN2 (Q(6) ,Q(7) 


TL ° DSIN(T3) 
T2 © DCOS(T3) 
pQq7) = 


T2*pQ(1) * 
T1*pQ(6) 


DQ (6) 
T2*DQ(6) - 
T1*DQ(1) 





exit 
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APPENDIX D 
Typical Computer Readout 

The following page represents a computer read-out of a particular 
transfer trajectory generated by the minimum time program previously 
described. It was generated at a step size of ten per revolution, i.e., 
step ten is the end of the first revolution. Read across from the step 
desired. The meaning of each column is indicated at its head. 

When reading across one will notice two lines of figures. In 
each case the second figure represents the change in value per radian of 
eccentric anomaly as the trajectory generates, with the exception that the 
second figure in the second column is the é function which must go through 
zero for optimality. On this page of sample read-out the ¢ function does not 
yet go through zero. But later on in the trajectory it begins to oscillate 


through zero with every revolution. 
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